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Abstract 

Radio Cherenkov radiation is arguably the most efficient mechanism for detecting 
showers from ultra-high energy particles of 1 PeV and above. Showers occuring in 
Antarctic ice should be detectable at distances up to 1 km. We report on electromag¬ 
netic shower development in ice using a GEANT Monte Carlo simulation. We have 
studied energy deposition by shower particles and determined shower parameters for 
several different media, finding agreement with published results where available. We 
also report on radio pulse emission from the charged particles in the shower, focusing 
on coherent emission at the Cherenkov angle. Previous work has focused on frequen¬ 
cies in the 100 MHz to 1 GHz range. Surprisingly, we find that the coherence regime 
extends up to tens of Ghz. This may have substantial impact on future radio-based 
neutrino detection experiments as well as any test beam experiment which seeks to 
measure coherent Cherenkov radiation from an electromagnetic shower. Our study is 
particularly important for the RICE experiment at the South Pole. 
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1 Introduction 

Ultra high energy (UHE) neutrinos can travel without scattering over large distances. 
These may prove to be useful cosmological and astrophysical probes. They also present 
themselves as candidate high energy particles with which we can test the Standard Model 
of electro-weak theory beyond the energy regime of current accelerators. In an UHE 
electron-neutrino charged current interaction, the neutrino gives most (~ 80%) of its 
energy to the secondary electron, which can then initiate an electromagnetic cascade or 
shower. It was predicted that an electromagnetic shower generated by a high energy 
primary could develop a charge excess which would emit Cherenkov radiation coherently 
0. For ultra high energy primaries the Cherenkov radiation would be coherent in the 
radio region of the spectrum; this long wavelength radiation might then be detected using 
radio antennas § I 1- Given the small predicted flux of ultra high energy neutrinos 
i !0> a suitable experiment to detect ultra high energy neutrinos using this method 
requires a large, dense (and radio-transparent) volume. Antarctic ice is suitable for this 
purpose. A detailed analysis of all aspects of such an experiment was done by Frichter, 


2 





































































Ralston and McKay (FRM) |§|, using a simulation developed by Zas, Halzen and Stanev 
(ZHS) FRM concluded that a modest array of optimally designed antennas could 
detect many events per year. The sensitivity of radio detection peaks above 1 PeV, which 
compliments an optical array such as AMANDA [p~0|l . 

The Radio Ice Cherenkov Experiment (RICE) at the South pole [11] is designed as a 
prototype detector of ultra high energy neutrinos with energy >PeV using this method. 
One basic requirement for such an experiment is a reliable Monte Carlo simulation of the 
shower development, Cherenkov radiation, detector, and data acquisition system. One 
can also test the idea of coherent Cherenkov emission at accelerator facilities by dumping 
a beam of photons or electrons into a dense target like sand or salt or any other suitable 
medium. Such tests have begun with experiments at Argonne and SLAC [12], A Monte 
Carlo simulation which can be easily adapted to such a test beam experiment, where 
Fresnel and possibly near zone radiation is important, and one that can include hadron 
showers conveniently, is clearly desirable. The ZHS simulation, designed for electromag¬ 
netic showers and Fraunhoffer (far zone) detection, has been a powerful tool. However 
an expansion and update with extensive testing, offering applications to test beam and 
neutrino astronomy, is currently needed. 

We have written a GEANT-based Monte Carlo simulation to study coherent Cherenkov 
emission in ice, salt, or a beam dump. GEANT 3.21 [|p|] is a well known and widely used 
simulation and detection Monte Carlo package in particle physics]*]. It allows access to 
all details of the simulation such as controls of various processes, definition of target and 
detector media, and a complete history of all events simulated. GEANT can be used to 
simulate all dominant processes in 10 keV - 10 TeV energy range, although it has not been 
extensively verified for energies >100 GeV, where the extrapolation of well-known lower- 
energy electromagnetic cross-sections becomes large, and other effects (LPM, e.g. 1B@) 
become significant. For electron energies above 10 GeV, GEANT uses screened Bethe- 
Heitler cross-section for bremsstrahlung together with the Migdal corrections [16, 17]. The 
first Migdal correction is important for energy >1 TeV, reducing the cross-section. The 
second correction reduces the differential cross-section for soft photon emission and is 
effective even at much lower energies, in the 100 MeV - 1 GeV range. The LPM effect 
in the context of UHE electromagnetic shower development and radio emission has been 
discussed in [[0], [li| ^0| . 

GEANT is used to simulate electromagnetic showers inside materials from which we 
extract detailed track information including particle type, coordinates, energy and inter¬ 
action time. From this track information, we investigate shower properties like radiation 
length, Moliere radius, critical energy and energy deposition in the material. We also de¬ 
termine the resulting radio pulse using standard electrodynamic calculations from charged 
particles’ tracks and by parametrizing the shower. Ultimately we will consider hadronic 
shower information and GEANT provides the flexibility to expand our analysis to this 
case. We note that physics results presented thus far by RICE have neglected the hadronic 
shower contribution and are, in this respect, conservative. 

differences between GEANT 3.21 and GEANT 4 are primarily at energies below the threshold for 
emission of Cherenkov radiation, and therefore do not affect the results presented here. 
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The organization of this report is as follows. In section two, we discuss various aspects 
of electromagnetic showers and define quantities which characterize the shower. We present 
results on the shower structure from the Monte Carlo simulation in section three and 
compare them with established values in standard materials such as iron, lead and carbon. 
Our analysis includes the detailed breakdown of the radiation-generating charge imbalance 
into energy bins, the direct evaluation from energy considerations of Moliere radius and 
the determination from dE/dx of radiation and energy deposition of the critical energy. 
We discuss shower-to-shower fluctuations by parametrizing the showers in section four. In 
section five, we review the theory behind coherent emission of an electromagnetic pulse 
from the shower [pi] ]. In section six, we calculate the electric pulse from the shower using 
track information from GEANT. We summarize our results in section seven, making a 
number of comparisons with ZHS, and discuss future work. 


2 Shower Description 

When a high energy electron or photon hits a material target, an electromagnetic shower 
develops longitudinally inside the material. At the beginning of the shower, bremsstrahlung 
and pair production are the dominant processes. Due to these processes, the number 
of particles increases exponentially and the shower is created. Due to the dominance 
of bremsstrahlung as the primary energy loss mechanism for high energy electrons, the 
population of photons quickly dominates that of electrons or positrons. The energy of the 
initial particle is divided among the secondaries. The exponential production of particles 
is halted when the charged particles reach the critical energy (E c ). This is the transition 
region where ionization loss overtakes radiation loss as the most important electron energy 
loss mechanism. The particle population also reaches its maximum (“shower max”) at this 
point. After reaching the critical energy, charged particles lose their energy predominantly 
inelastically by ionization, resulting in a subsequent decline in the number of particles in 
the shower as they degrade in energy and range out in the medium. 

Multiple Coulomb scattering is responsible for the transverse spread of the shower. The 
shower core is populated by the highest energy particles. There is a long tail composed 
mostly of the coulomb- scattered low-energy component. 

Other processes - Compton, Moller and Bhabha and positron annihilation build up a 
net charge (more electrons than positrons) in the shower as atomic electrons in the target 
medium are swept up into the forward moving shower. 


2.1 Energy Loss Mechanisms 


Energy loss of an electron due to radiation is well approximated by the Bethe-Heitler 
formula @ 



= E 


m o Z 2 rj 
137 A 


In 



(1) 


where E is the initial energy of the electron, x is the distance in g-cm 2 units, N 0 is 
Avogadro’s number, A is the mass number of the nuclei of the medium, Z is the atomic 
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number of the medium and r 2 = e 4 /m 2 c 4 is the classical electron radius. For a composite 
medium like ice (HyO), one has to calculate the effective Z and A using the proportion 
by weight method 2 ]. 

] gives the energy loss due to ionization as 


The Bethe-Bloch formula [ 23 


d E\ 

) ion 


-r- I = AirN 0 —r/mc 2 


ln /2m^V\ _ 1 _ 6 


( 2 ) 


where I = 10 Z eV and 7 is the Lorentz factor in usual relativistic notation. At very high 
energy, one needs to take into account the density effect as the medium becomes polarized. 
This is taken care of by the following term inside the square bracket of Eq. (I): 


5/2 = \n(hu) p /I) + ln/Ty — 1 / 2 , 


( 3 ) 


where w p is the plasma frequency of the medium. 

A rough expression for the critical energy, E c , can be found by the ratio of the expres¬ 
sions for radiation and ionization losses. The log term in Eq. (jl|) is roughly 4 and the 
square bracket term in equation (j]|) is approximately 11 if we include polarization effects. 
The ratio is then 


E c 


605 

IT 


MeV. 


( 4 ) 


An alternate definition of the critical energy is that of Rossi [ 251, according to which, 
the critical energy (E c ) is the energy at which the ionization loss per radiation length is 
equal to the electron energy. Critical energies determined from these two definitions are 
compared in a subsequent section. 

The radiation length (X 0 ) is given by 


1 4 N 0 Z 2 r 2 e / 183 \ 

X 0 ~ 137A n V^ 1/3 y 


( 5 ) 


so that the approximate expression for radiation loss (Eq. ®) becomes — (dE/dx) = 
E/X 0 or 

<E>= E 0 e~ x/X °. ( 6 ) 

Eq. (j 6 [) serves as a definition of the radiation length; X a is the thickness of material needed 
to reduce the mean energy of an electron to 1 /e of its original through bremsstrahlung. 


2.2 Longitudinal Profile 

Heitler |2(| developed a simplified model of an electromagnetic shower according to which 
the initial electron, with energy E 0 , radiates a photon, of energy E Q / 2, in the 1st radiation 
length. In the 2nd radiation length, the photon creates an electron-positron pair and 
the previous electron emits another photon; each particle now has energy E 0 / 4. This 
process continues until the critical energy E c , defined in the previous section, is reached, 
at which point the shower is defined to be at maximum. At any radiation length, t, the 

2 see for example CONS110 section of Jl3| . 
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number of particles (electron, positron and photon) is N(t) = 2 t = exp (tin2) and the 
energy per particle is E(t ) = E 0 /2 t . Thus in this model the shower maximum occurs 
at tmax = In (E 0 /E c )/\n2 and the total number of particles at the maximum is tV max = 
E 0 /E c . We can also calculate the integrated track length of the charged particles as 
L = (2/3ln2)(E 0 /E c ). 

Although this model is over-simplified in assuming equality of all the primary cross- 
sections and in assuming that the shower is artificially cut off, it nevertheless provides 
a good description of the qualitative features of an electromagnetic cascade. Namely, it 
predicts that the total track length of the charged particles scales linearly, while the posi¬ 
tion of the shower maximum scales logarithmically with the initial energy of the shower. 
These are important consistency checks for the electromagnetic showers generated from 
any Monte Carlo. 

An analytic method for realistic shower development was created by Carlson and Op- 
penheimer |27j and later extended by Landau and Rumer |28|] . Analytic shower equations 
are solved in two approximations, namely approximation A and approximation B 
One neglects the Compton effect and ionization loss in approximation A. One includes a 
constant energy loss per radiation length for all electrons and positrons in approximation 
B, which is therefore expected to be a better model of the data. The derivation of the 
integral spectra or number of particles (electrons and positrons) at different shower depthf^j 


t in approximation A and B can be found in references [25, p9| . 

Greisen first parametrized the longitudinal profile of a photon-induced electromagnetic 
air shower |h|. The Greisen parametrization (GP) is based on analytic shower theory in 
approximation A. The difference between the more realistic approximation B and approx¬ 
imation A is that of a slowly increasing function of the age parameter (s). The number of 
particles in a given energy range increases when s < 1 , reaches a maximum when s = 1 
and then declines when s > 1. 

Hillas [|TJ later modified Greisen’s parametrization to fit Monte Carlo simulation of 
0.1-1 TeV photon-induced showers. This modified GP was later used by Fenyves et. al. 
and is given for a single electron or photon of energy E 0 as: 


N^\E a ,E,t) = exp[fi(l - 1.5 In si)] 

y/y 


( 7 ) 


where t\ is the modified depth defined to be t\ = t + a v ^{E) and si is the modified age 
parameter defined to be si = 2>ti/(t\ + 2 y). Superscript n denotes the total electrons and 
positrons following the convention in the literature. The parameter y called lethargy is 
defined to be y = In (E 0 /E c ), where E c is the critical energy. One finds the parameters 
A(E), a n (E) for electron-induced and a-y(E) for photon-induced showers by fitting Monte 
Carlo simulations. 


2.3 Lateral Spread 

The transverse development of the shower is well described by a quantity called Moliere 
radius (Rm)- It is determined by the average angular deflection per radiation length at the 

3 Shower depth t = depth/N 0 is always given in terms of the radiation length (X Q ). 
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critical energy (E c ) due to multiple Coulomb scattering. The average deflection is given 

<^ 2> =(§) 2 ^ ( 8 ) 

where E s = \J47 r /a m e c 2 = 21.21 MeV is the scale energy. 

A numerical estimate of Rm is obtained by considering the fraction of energy that 
escapes transverse to the shower axis [|34], [ll| : 

U(r) = /q° / r °° E(r, t)drdt 
E /o°° fo°° E ( r . t)drdt 


where energy (E) is expressed as a function of shower depth (t) and radial distance r from 
the shower axis. By definition, ninety percent of the shower energy is contained inside a 
cylinder of radius Rm centered on the shower axis. I.e., r = Rm when U{r)/E = 0.1 in 
Eq. Moliere radius is independent of the energy of the shower and depends only on 
the tracking medium in general. 

In Rossi’s definitions, the Moliere radius is related to the critical energy (E c ) and 
radiation length (X 0 ) of the material [[k|, |35|] through the equation 


Rm 


XqEs 

Ec 


( 10 ) 


which follows from Eq. (||). 


3 Shower Simulations 

The target medium in our GEANT simulations, unless otherwise stated, is defined to be 
an ice cube of side 1 kilometer. Given the molecular composition, GEANT calculates the 
effective atomic number, Z = 7.2 and an effective mass A = 14.3 for the compound ice. 
Other parameters like radiation length, absorption length and cross-sections are calculated 
automatically, once A and Z are specified. 

GEANT gives all the details of particle tracking information like interaction points, 
total energy, energy lost in interaction, interaction time and so on. We used double 
precision for all the variables in our output data files to minimize roundoff errors. Unless 
stated otherwise, all particles are tracked down to total energy of 0.611 MeV, which is lower 
than the energy at which particles are still relativistic and emit Cherenkov radiation. 

3.1 Radiation Length 

The radiation length given by Eq. (||) depends on the atomic and mass numbers of the 
material. However employing the definition given by Eq. (|[), we can fit an exponential to 
the average energy of the injected electron at increasing depths. The average is taken from 
the data of a number of showers generated by the Monte Carlo. Our ability to recover the 
input value of X D will serve as an internal consistency check to ensure that we are tracking 
all the particles of interest, along with their energies. 
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Figure 1: Average energy of a 1 TeV electron injected in ice vs. distance in cm. Monte Carlo 
data points are obtained from 500 electron tracks and the errorbars correspond to standard error 
or s/y/N, where s is the standard deviation and N is the number of tracks (500). The solid line 
is the exponential fit as in Eq. (D with E a = 1 TeV. This fit gives a radiation length of 42.2 ± 4.3 
cm. A 100 MeV energy threshold was used in all 500 simulations. 


We generated 500 electron tracks each with 1 TeV primary energy ( E rmo ) incident 
on the ice target and calculated energy loss due to bremsstrahlung. We tracked all the 
electrons down to energy 100 MeV, well above the value of the critical energy so that 
bremsstrahlung is dominant. We show the energy of the injected electron, averaged over 
the 500 tracks, as a function of distance in Fig. [I]. The errorbars correspond to standard 
error or s/y/N , where s is the standard deviation and N is the number of tracks (500). 
A least squares fit of Eq. (I) to the Monte Carlo data keeping fixed E 0 = 1 TeV yields a 
radiation length of 42.2 ± 4.3 cm. The confidence level (CL) for the fit is 95.7%. 

Given the molecular composition, GEANT also calculates the medium’s radiation 
length from the standard formula. For ice, GEANT calculates X Q = 38.8 cm, which 
is roughly in agreement with the value we extract by tracking bremsstrahlung photons. 

3.2 Moliere radius 

To calculate the Moliere radius, we construct an imaginary cylinder centered on the shower 
axis (up to the physical length of the shower). We add the energies, U, of all the tracks 
that leave the cylinder without re-entering. By varying the radius of the cylinder we obtain 
Rm, which is the radius of the cylinder when the fraction U / E Q is equal to 0.1. 

The Moliere radii (Rm.) f° r lead and iron are found to be 1.6 cm and 2.1 cm, respectively 
(see Fig. ^). The value for lead is in agreement with experimental value of 1.6 cm [|34|]. 
We did not find the corresponding value for iron in the literature. A similar analysis for 
ice results in R^ = 13 ± 1 cm. We quote the standard error, defined as s/y/N, where s is 
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Figure 2: Moliere radius corresponds to the transverse development of the shower. When the 
fraction U/E a (the ratio of total energy inside an imaginary tube along the showers axis to the 
initial energy of the shower) is 0.1 (the horizontal straight line), the corresponding radius is by 
definition, the Moliere radius for that material. It depends on the material and not on the energy 
of the shower. Here we compare showers of energies 1 GeV (left plot) and 10 GeV (right plot) for 
each material; the Moliere radius is observed to be independent of energy in this range. 


the standard deviation and N is the number of showers. It should be noted that the error 
bars are correlated from bin to bin. 

3.3 Energy Loss in the Medium 

The signal from a Cherenkov type detector is proportional to the track length, which 
is itself proportional to the energy deposition of the shower particles in the medium. 
Therefore, we study the ionization loss per unit length in this section. 

To determine dE/dx due to ionization from Monte Carlo, we generated 500 separate 
5 GeV electron tracks in ice and kept a record of the rate at which energy was lost due 
to ionization. The dots in Fig. |3] (left plot) shows dE/dx in ice for the 500 tracks at 
different Monte Carlo steps. We also calculated dE/dx ionization loss in carbon, which 
has roughly the same atomic weight as ice, as a consistency check. The result is shown in 
Fig. H (right plot). 

The average dE/dx curve from GEANT simulation matches the analytic approxima¬ 
tion, Eq. (p) as shown in Fig. [|. The average ionization loss in ice (left plot) in the 
relativistic region is approximately 2.4 MeVg -1 cm -2 and the average minimum ioniza¬ 
tion loss is approximately 1.9 MeVg -1 cm -2 . We have also calculated the average dE/dx 
ionization loss in carbon (right plot) as a consistency check. 
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Figure 3: Electron stopping power or dE/dx ionization loss in ice (left plot) and in carbon (right 
plot). The dots correspond to energy lost per unit length at each Monte Carlo step of 500 electron 
tracks each with 5 GeV primary energy. 
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Figure 4: Average electron stopping power or dE/dx ionization loss in ice (left plot) and in 
carbon (right plot) calculated from GEANT data in Fig. ||. The error bars correspond to standard 
error or s/y/N, where s is the standard deviation and N is the number of tracks (500 in this case). 
The solid lines are analytic curves calculated from Eq. (|j). The agreement between Monte Carlo 
and the analytic formula is reasonably good in the relativistic rise region. 
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3.4 Critical Energy 


There are several methods to calculate the critical energy (E c ). The simplest is of course, 
the rough estimate in Eq. ([|). One can also obtain a value for E c in the Rossi definition 
from Eq. & after supplying the value of the Moliere radius (Rm) and the radiation 
length (X 0 ) found earlier. Yet another method is to find the energy at which the dE/dx 
curves due to ionization and radiation losses intersect, which is the Monte Carlo equivalent 
of the analytic expression in Eq. (||). 

First, we find the critical energy from the Rossi definition using Eq. ©■ The Moliere 
radius for ice was found to be = 13 ± 1 cm and the radiation length, X Q = 42.2 ± 4.3 
cm. The critical energy is then calculated from Eq. (10) to be E (f 6 = 68.8 ± 8.8 MeV. 
Rossi [^] quotes 65 MeV for the critical energy of water, in agreement with our calculated 
value. As a consistency check, we also calculated the critical energy for lead using the 
Moliere radii found from GEANT simulation (section 3.2). The result is E^ ead = 7.4 MeV, 
where we have used the radiation length of lead to be 0.56 cm. This is in agreement with 
the experimental result [36|. 

Second, we calculate the energy loss due to radiation according to the formula in Eq. 
(§) using the radiation length X 0 = 42.2 ± 4.3 cm we found before. We then find the 
energy at which this line crosses the ionization loss curve (see Fig. |^). The value for 
critical energy (the energy at the crossing point of these two curves) found in this method 
is 90 ± 9 MeV, where the error includes both the standard error in the dE/dx points and 
the uncertainty in the radiation loss due to the error in X Q . The approximate formula 
(Eq. (||)) gives E c = 84 MeV, close to what we found from simulation. The critical energy 
calculated in this way gives a higher value than Rossi definition, as expected. 


3.5 Track Lengths 

If we consider all processes to be elastic except ionization, then an estimate of the upper 
bound for the total track length is given by, 


L = EJ 


(dE 

\d x J ■, 


( 11 ) 


where (dE/dx)™™ is the minimum ionization energy loss per unit length. This value is 
about 1.9 MeVg -1 cm~ 2 for ice, (see Fig. [|) which yields a maximum total track length 
of ~ 1350 radiation lengths or 570 meters in ice for a 100 GeV shower. The actual total 
track length is less than this value, since the energy loss due to ionization increases in the 
relativistic region. 

We plot different track lengths in Fig. [6] for a 100 GeV shower (averaged over 25 
showers) versus the kinetic energy threshold used in the Monte Carlo to generate them. 
The total track length for electrons and positrons together is denoted by total absolute 
track length, the sum of electron and positron track-lengths projected along the shower 
axis is denoted by total projected (e + p) track length and the difference between the 
electron and positron track lengths projected along the shower axis is denoted by total 
projected (e — p) track length. As we increase the threshold, so the simulation neglects 
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Figure 5: Radiation and ionization energy loss vs. electron energy plot from GEANT Monte 
Carlo simulation. The diagonal straight lines correspond to energy loss due to radiation with 
radiation length X Q = 42.2 ± 4.3 cm, where the lowest curve corresponds to X Q = 42.2 + 4.3 cm 
and the highest to X a = 42.2 — 4.3 cm. The average ionization energy loss data are plotted as 
obtained before (see Fig. ||). The error bars show the standard errors on the dE/dx points. The 
crossing point of the two curves gives one definition of the critical energy. 


track-lengths from low energy particles, we expect these total calculated track-lengths to 
decrease. 

We have also compared our results to those obtained using the ZHS Monte Carlo in 
Fig. |j. Although the qualitative behavior of this scaling is the same for both Monte 
Carlos, ZHS track lengths are consistently higher than GEANT track lengths by about 
50%. In particular, GEANT produces an absolute track length of 400 meter for 0.1 MeV 
kinetic energy threshold or 0.611 MeV total energy threshold. ZHS on the other hand 
produces 650 meters of absolute track length for the same threshold, which is significantly 
above our estimated maximum of 570 meters, based on the GEANT generated 100 GeV 
data shown in Fig. || and application of Eq. (0. 

The total track lengths increase with shower energy as more particles are created. The 
track lengths are expected to increase linearly with energy for a given a energy threshold 
used in the Monte Carlo. This scaling is shown in Fig. [i] below. Straight line fits to the 
absolute track length and the projected (e — p) track length are also plotted. The slopes 
of those straight lines are 3.2 and 0.5 respectively. 

3.6 Particle yield 

To generate a longitudinal profile of the shower, one counts the number of particles crossing 
planes perpendicular to the shower axis inside the ice. We have calculated profiles for both 
the total number of particles (e + p) and for the excess charge (e — p) in the shower. 
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Figure 6: Total track lengths as functions of the kinetic energy threshold used in the Monte 
Carlos. Here we show the total absolute track-lengths, the total (e + p) track-lengths projected 
along the shower axis and the total (e — p) track-lengths projected along the shower axis. The 
analysis is done for a 100 GeV shower (averaged over 25 showers) using the GEANT and ZHS 
Monte Carlo simulations. 
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Figure 7: Total absolute track length, total projected (e + p) track length and total projected 
(e —p) track length as functions of shower energy from GEANT simulations. The linear scaling of 
the track lengths with the shower energy is clearly observed. We used a kinetic energy threshold 
of 1 MeV (see Fig. ^) and averaged over 25 showers in each case. 
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Figure 8: The longitudinal profiles of a 30 GeV shower (averaged over 50 showers) in iron 
simulated by GEANT. The longitudinal profile is calculated by counting the number of electrons 
and positrons crossing planes on the shower axis every half radiation length with total energy 
greater than 1.5 MeV. We also plot the longitudinal profiles of the same shower from EGS4 (from 
PDG 2000 J37t ) and by ZHS Monte Carlos. 


As a consistency check of GEANT, we have simulated 30 GeV electron-induced showers 
in iron and compared the profiles with the modified Greisen parametrization, Eq. 

The longitudinal profiles were obtained by adding the number of particles with total 
energy greater than 1.5 MeV, crossing planes spaced one-half radiation length apart, and 
perpendicular to the shower axis (see Fig. ^). The number of particles agrees reasonably 
well with EGS4, an electromagnetic shower code developed at SLAC which simulates the 
same shower [37]. The total number of particles from ZHS simulation is also shown in the 
plot. 


We have fitted the GEANT generated electromagnetic showers in ice by the modified 
GP in Eq. ©■ We used the value for critical energy ( E c ) in ice to be 67.7 MeV. The 
fitting parameters and the confidence level (CL) of the fits for showers of energy 100 GeV 
(averaged over 100 showers), 500 GeV (averaged over 50 showers) and 1 TeV (averaged 
over 20 showers) each with 5 MeV threshold energy are given in Table 1 for both electron 
and photon primaries. The fits to longitudinal profiles are plotted in Fig. |9|. 
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Figure 9: Longitudinal profiles of electron-induced (left) and photon-induced (right) showers 
in ice of primary energy 1 TeV, 500 GeV and 100 GeV, averaged over 20, 50 and 100 showers 
respectively. The threshold energy used is 5 MeV. The error-bars correspond to standard error. 
The smooth curves are Greisen parametrization fits as given in Eq. ([?]). The fitting parameters 
and the confidence levels (CL) for the fits are given in Table 1. The critical energy for ice used in 
the fits is 68.8 MeV. 


Table 1: Greisen parameters for GEANT 
showers in ice with 5 MeV threshold and 
E c = 68.8 MeV. 


Primary 

E 0 (GeV) 

A{E) 

^7r, 7 

CL(%) 


100 

0.50 

0.33 

96.3 

7 

500 

0.50 

0.26 

81.9 


1000 

0.50 

0.76 

71.8 


100 

0.52 

0.99 

95.6 

e - 

500 

0.51 

1.14 

95.8 


1000 

0.52 

1.01 

88.4 


The parameter A(E) gives over-all normalization for the number of particles. We see 
from Table 1 that the number of particles at the shower maximum for a photon-induced 
shower ( A(E ) = 0.50) is slightly less than that of an electron-induced shower (A(E) = 
0.52). The other parameter takes care of the shift of the whole shower between a 
photon-induced (a 7 ) and an electron-induced (a^) shower. This shift is approximately 0.7 
radiation length. 

We have also fitted GEANT generated 100 GeV electron and photon induced showers 
(each averaged over 50 showers) in ice with 0.611 MeV threshold energy by GP given in 
Eq. ©■ The fitting parameters and the CL for the fits are given in Table 2. The fits with 
longitudinal profiles are plotted in Fig. 0 The critical energy in ice for the fits is 68.8 
MeV as before. 
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Figure 10: Greisen parametrization fits (smooth curves) to the longitudinal profiles of 100 GeV 
electron and photon induced showers (each averaged over 50 showers) with 0.611 MeV threshold 
energy from GEANT. The fitting parameters are given in Table 2. The critical energy for ice used 
in the fits is 68.8 MeV. 


Table 2: Greisen parameters for GEANT 
showers in ice with 0.611 MeV threshold 
and E c = 68.8 MeV. 


Primary 

E 0 (GeV) 

ME) 

C^7 r, 7 

CL(%) 

7 

100 

0.65 

0.16 

85.6 

e~ 

100 

0.66 

0.66 

89.1 


A three parameter fit to the GEANT generated 100 GeV electron-induced shower in 
ice leaving the critical energy (E c ) as a free parameter along with A{E) and a n gives 
E c = 58.3 ±3.1 MeV. This is within error bars of the critical energy (E c = 68.8 ± 8.8 
MeV) found from Moliere radius calculation. The CL for this fit is 97.7%. 

A comparison of averages over 50 showers of 100 GeV each from GEANT to the same 
from ZHS Monte Carlo (see Fig. 0 ) shows about a 25 — 35% discrepancy for the total 
number of particles at the shower max. The percentage discrepancy between the two 
simulations remains the same at higher energies. 

The excess charge in the shower (N(e — p)/N(e + p)) is about 18% at the shower 
maximum, about twice Askaryan’s original rough estimate JTj. 

Fig. H shows a longitudinal depth distribution of a much different kind. Unlike 
particles crossing planes transverse to the shower axis, here we look at the distribution 
of the complete count of excess electrons present in bins of one radiation length along the 
shower axis. The figure shows the distribution of excess electrons (dlV(e — p) /d t) in a 
100 GeV shower (averaged over 50 showers) with 0.611 MeV total energy threshold. The 
distribution is broken down into different energy bins, showing that half of the excess 
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Figure 11: Comparison between shower profiles of a 100 GeV shower (averaged over 50 showers 
in each case) from GEANT and from ZHS Monte Carlos. The error bars correspond to standard 
error, or s/y/N, where s is the standard deviation and N is the number of showers. The total 
energy threshold used is 0.611 MeV in both cases. The difference in particle yield between the two 
Monte Carlos is 25 — 35% and remains the same at higher energy. 


electrons have energy 5 MeV or lower. This plot shows clearly the important point that 
the bulk of the particles at shower maximum have low energy and contribute significantly 
to the Cherenkov emission. 


3.7 1-D Model Study 


The comparisons in Figs. |g & 11 prompt us to check the sensitivity of shower depth and 
particle number at maximum to the cross-sections assumed in a simulation, since ZHS 
and GEANT codes differ slightly in cross-sections for some processes (a few percent). To 
investigate this question in a setting that is completely under our control, we have written a 
one dimensional Monte Carlo that elaborates on Heitler’s simple model to make it realistic. 
We include Compton scattering, positron annihilation, ionization loss, Moeller scattering, 
where interaction lengths are chosen from exponential distributions. The default values 
for the cross-sections are taken from GEANT. We examine the effect of changing cross- 
sections of various processes. The total number of particles (e + p and 7 ) from the 1-D 
model scales linearly with shower energy and the position of the shower maximum scales 
logarithmically with energy as expected. 

The effect of changing individual cross-sections (increasing and decreasing by 25%) are 
listed in Table 3 & 4. We also list the total number of electrons and positrons (N(e + p)) 
at the maximum depth (tffififi) and the total number of photons (V(y)) at the maximum 
depth (t'h ax ) for the default values of the cross-sections (denoted by a Q ). 


17 









Depth (t) 

Figure 12: Depth distribution of excess charge (d N(e — p)/dt) broken down into 3 different 
energy bins: E < 5 MeV, 5 < E < 50 MeV and E > 50 MeV. This figure shows that half of the 
excess electrons are created with energy below 5 MeV. We used a 100 GeV shower (averaged over 
50 showers) with 0.611 MeV threshold to generate this plot. 


Table 3: Results of increasing individual cross-section by 25% in the 
1-D Monte Carlo. The default setting (GEANT) is denoted by a D . 

The cross-sections are changed one at a time from their default values. 


Parameter 

cr 0 

Brem. 

Pair 

Comp. 

5 

dE/dx 

Anni. 

f e+p 

‘'max 

6.9 

7.2 

6.2 

6.9 

5.6 

6.4 

6.7 

N(e + p) 

266 

236 

307 

281 

269 

231 

265 

f1 

"max 

6.9 

7.0 

5.9 

6.4 

6.4 

6.7 

6.4 

N( 7) 

1418 

1570 

1495 

1292 

1396 

1209 

1416 


Table 4: Results of decreasing individual cross-section by 25% in the 
1-D Monte Carlo. The default setting (GEANT) is denoted by a Q . 

The cross-sections are changed one at a time from their default values. 


Parameter 


Brem. 

Pair 

Comp. 

(5 

dE/dx 

Anni. 

j.e+p 

"max 

6.9 

6.4 

7.7 

7.4 

7.4 

7.4 

6.1 

N(e + p) 

266 

310 

212 

249 

298 

298 

271 

fl 

"max 

6.9 

6.7 

8.2 

7.4 

6.9 

6.9 

6.4 

N( 7) 

1418 

1173 

1259 

1516 

1544 

1544 

1399 


The individual cross-sections are changed one at a time from their default values. The 
changes in the values of t max and N are noticeable, but the percentage changes are much 
smaller than the percentage changes in the cross-sections. For example, comparing the ef¬ 
fect of increasing the bremsstrahlung cross section, listed in the column under “Brem.” in 
Table 3, the depth and number at maximum change by at most 4% and 15% respectively, 
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compared to the default values, listed under “<r 0 , for a 25% change in the cross section. 
It seems that the shower features are not especially sensitive to precise values of an in¬ 
dividual cross section. We have not made an exhaustive study of the effects of changing 
combinations of cross sections in all possible ways, but we have seen no indication that 
effects would be large without pathological distortions of the cross-sections. 


4 Shower Fluctuation 


Any determination of shower energy is intrinsically uncertain because of the fluctuations 
from shower to shower for the same injected energy. The fluctuations must be quantified 
compactly so an efficient and realistic energy uncertainty from this source can be assigned 
to detected showers. Fortunately, for UHE showers the huge numbers of particles and 
fractionally smaller fluctuations make this less problematic than for small showers. We 
summarize an approach to the fluctuation question in this section. 

Conservation of total energy plays an important role in shower development. The shape 
of a shower strongly depends on the position and energy of the first hard bremsstrahlung - 
the later the first hard bremsstrahlung event, the deeper the shower maximum is. However, 
for any shower, the primary energy (i.e., the total energy of the shower) qualitatively 
dictates the energy loss profile and the particle yield with depth. 

It is common practice to describe the mean longitudinal profile of an electromagnetic 
shower by a Gamma distribution which is similar in shape to the Greisen parametrization. 
The Gamma distribution normalized to unity is given by 




b 


( bt) a 1 e 



( 12 ) 


However, Grindhammer et. al. (3^] have shown that the fluctuations of the parameters 
a and b from an average profile do not necessarily follow the individual shower fluctuations. 
It is more reasonable to fit individual profiles by the Gamma distribution. One can then 
fit a 2-dinrensional Gaussian distribution to the parameter set {a, b} thus obtained and 
study shower fluctuation by studying these parameters. The correlation between these 
parameters a and b can be expressed as: 


p = Covariance (a, b)/a a crb, 


(13) 


where cr a and are standard deviations of the parameters a and b respectively and 
Covariance (a, b ) is the covariance matrix. The correlation p is roughly independent of the 
energy of the shower. 


4.1 Fit to Simulations 


We have fitted the Gamma distribution Eq. (12) to 50 GEANT generated individual 
shower profiles in iccQ. The energy of each shower is 100 GeV with 0.611 MeV threshold 
energy. The parameters a and b are extracted from the Gamma function fits. Two such 


fits to individual shower profiles are plotted in Fig. 13 
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Shower Depth Shower Depth 

Figure 13: Individual shower profiles of 100 GeV showers from GEANT with 0.611 MeV threshold 
energy. The solid curves are Gamma function fits to the profiles. The parameters a and b in Eq. 
( 0 ) are obtained from the fits. 




Figure 14: Scatter plot (left plot) of the parameters a and b from the Gamma distribution fit to 
50 individual shower profiles each with energy 100 GeV and 0.611 MeV threshold. The contours 
are from a 2-dimensional Gaussian distribution with mean (a, b) and (cr 0 , cq,) obtained from the 
data set a, b generated as described in the text. Shower fluctuations (right plot) due to variation 
of the parameters a and b within a standard deviation. The dark solid curve is the profile with 
mean values of a and b. All particle numbers are normalized to 1. 
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Fig. |l4j a shows the scatter plot of the parameters a and b obtained by fitting 50 
individual showers. We also show the contours of the 2-dimensional Gaussian distribution 
fit to the set of parameters {a, b}. This plot is similar in shape to the one obtained by 
Grindhammer using GEANT |3£|. The correlation coefficient p = 0.72 we found is close 
to 0.73 as found by Grindhammer et. al. 

Fig. 0b shows shower-to-shower fluctuations as we vary the parameters a and b within 
a standard deviation ( a a and cp,) of their mean values ( <a> and <b>). The dark solid 
curve corresponds to the shower with mean values of the parameters a and b. As can be 
seen, the shower maximum can vary as much as ~ 50%. As remarked above, the variations 
become statistically less significant as the energy rises to the multi-PeV range, where UHE 
neutrino detectors such as RICE are sensitive. The uncertainties from fluctuations, though 
requiring study, are not expected to be a major part of the uncertainties in event locations 
and energies in any case. 


5 Electromagnetic Pulse Theory 

Overview 

Given a current density, one calculates the electric field by straightforward application 
of Maxwell’s equations. The most useful form of the electric field is the complex vector 
E w (x), depending on three-dimensional coordinate x and angular frequency oj. Before 
getting into the details of the radiation field calculation, we review its general features 
and several subtle points involving coherence, and near and far field limits. 

The energy radiated per unit frequency interval per unit solid angle is proportional to 
l-E^I 2 . In an asymptotically far field limit R —> oo, the radiation field E w oc 1 /R, where R is 
the distance from the shower to the field point. Because the quantity E u has dimensions 
of mass in particle physics units and depends linearly on the track length L (reviewed 
below), one expects the frequency dependence E u oc Lu/R on dimensional grounds. This 
linear dependence on c o breaks down at high frequencies, when wavelengths are smaller 
than the typical dimensions of the shower and create a coherence cutoff. The purpose of 
this section is to develop the concepts and formalism needed to address this and other 
issues. 

The asymptotically far field is not a feature of Cherenkov physics as commonly pre¬ 
sented in texts or applied in particle physics detectors. In particle physics applications, 
the track length is comparable to the distance to the observation point R. Textbook 
treatments take the limit that the track length is infinite, [39, 4C, in which case the 
electric field E w ~ y/uj/^/~p, where p is the cylindrical radial coordinate. In this situation, 
the radiation intensity is proportional tow, a familar feature of laboratory Cherenkov ra¬ 
diation. The cylindrical symmetry dictates this special frequency and radial dependence. 
The transition from cylindrical configuration to spherical configuration is associated with 
the terms “Fresnel and Fraunhoffer” zones. (Conditions of “far fields” R A, p A 
are separate and assumed throughout; criteria separating the zones are given below.) The 


4 particle number (e + p) normalized to 1. 
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conditions of the RICE experiment are such that spherical symmetry is a good approxi¬ 
mation for the most distant events, R ~ 1 km, while Fresnel zone effects start to become 
noticeable for the nearest likely events, R ~ 100 m. Consequently both are reviewed. 

There are two methods to calculate the electric field from the charged particles in a 
shower. One method calculates the electric field from each charged track and adds them 
by superposition. The other method parametrizes the shower’s current and calculates elec¬ 
tric fields analytically. Direct track-by-track calculation and superposition are presently 
limited to the Fraunhoffer zone. The result is simple, and depends only on the angle to 
the observation point. The analytic method works in either zone, and gives a compact 
way to take a few parameters characterizing the shower and calculate the field. Of course 
one must have a good parameterization for this to be a good approximation. 

By using both the analytic and numerical methods, with their different strengths, we 
are able to characterize the coherence structure of showers in considerably more detail 
than previous work. The coherence extends to much higher frequencies than previously 
thought. This is a very important point, confirmed by two independent methods, and 
clears up a misconception that coherence should cut off at about 1 GHz. 


5.1 Fraunhoffer Limit 


The Fraunhoffer limit is appropriate to most of the RICE sensitive volume. It forms the 
basis for our numerical, track-by-track computation of the field produced by a shower at 
antenna sites remote from the shower. In this subsection, we derive the expression for 
the field produced by an individual track, which forms the basis for the calculation of the 
full field calculated from a whole shower. In the next subsection, we outline the method 
of parametrizing the shower as an effective current and calculating the field directly from 
that current. 

The power at time t, position x, radiated per unit solid angle by a moving electric 
charge is given, in Gaussian units, by |39fl 

”” d?P = -^\ R E( x it)ret\ , ( 14 ) 

evaluated at the retarded time t = t' + R{t')/c. (For now we calculate effects in vac¬ 
uum; shortly we will supply the factors for the effects of a medium.) We use the Fourier 
transformed variables 

] roc 

REu = —= / RE(t)e lwt dt ; (15) 

v 27 T J — oo 

_ 1 roc 

RE(t) = —= / RE u e~ lut duj . ( 16 ) 

y/2ir J —oo 

The energy radiated per unit frequency interval per solid angle is then given by 


d 2 I 

dtudG 



(17) 


We will require the frequency dependence of the fields, so we work with the Fourier trans¬ 
formed fields below. 
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Figure 15: Geometry for calculating electromagnetic fields from a single track segment. (ri,fi) 
and ( 72 ,^ 2 ) are the starting and ending position and time of the segment along which the particle 
moves with velocity v. 


The expression for the radiation field from a point source is conventionally defined as 
the electric field term linear in the acceleration (3: 


E(x, t) = - 
c 


h x {(n — (3) x /3} 
(1 -/3-h) 3 i? 


(18) 


ret 


where (3 is the velocity of the particle, h is the direction of the observer and R is the 
distance from the track to the observation point (see Fig. |l5|) . The factor 1/R that 
accompanies the (3 factor is the other trademark of the radiation held. As is the case for 
the term that comes from the boosted Coulomb held, which has no explicit acceleration 


dependence, Eq. (18) is singular at the Cherenkov angle in a medium with real index of 
refraction greater than 1. 

Combining Eq. fll5l ) and Eq. ([l^) we have 


Eujix) = 


/ r Mt+R(t)/c) 

fi x {(h — (3) x (3} 

/ 87T 2 C J- 00 

1 

1 

3 

to 

ba 

l_ 


d t'. 


(19) 


At distances large compared to the range of motion of the source, h is approximately 
constant and 


R(t') « |a?| — n • r{t'). (20) 

This is the Fraunhoffer approximation : the error in the phase oj\x — x'\/c must be kept 
small compared to 27r. Conditions for use of this approximation are discussed in the next 
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subsection[] (see Eq. (|29j) ). 

After integrating Eq. (H) by parts, and using the boundary conditions to set the end 
point contributions to zero, one finds [39| 


RE U (x) « -ieo\ ^e iu}R/c r e Mt-n-r/c)^ x ^ x 
V 07T^C 7-oo 


( 21 ) 


where i? = |x|. We will later apply this track-by-track expession to segments over which 
fi is constant with time from t\ < t < t 2 - On each segment r = fq + cfi(t — t\). We insert 
the result into Eq. (21) and obtain 


RE u] (x) = 




^ e ikR e iui(ti-h-ri/c)p 


^ e iwSt(l—h-0) _ 

1 — n • fi 


( 22 ) 


where h x (n x fi) = — fix, and St = t 2 — 1\. 

In a medium we replace c —> c/n = c/yfefi. where n = n( u) is the refractive index 


of the material and e = e(o;), the dielectric constant. We also replace E 
(3 —► = nf3 everywhere, to get 


RE^(x) = 


1 


_\ ( c iui8t(l—h-f3n) 


1 ) 


\Z2tt \ c 2 J 


1 — n ■ /3n 


D = eE and 


(23) 


where /j, r is the relative permeability and k = nuj/c. 

The particle velocity v in the medium can be greater than that of light. The apparent 
singularity in Eq. (|2^) at 1 — n ■ fin = 0 defines the Cherenkov angle Q c as cos 9 C = 1/nfi. 
However there is no singularity, as seen by expanding the exponent. Close to the Cherenkov 
angle, Eq. (|23|) reduces to the form 


RE u (x) = -^= (^f) e^e^-^'^vxSt. 
V27T V C Z J 


(24) 


Eqs. (E3l&;l24|), used by ZHS without detailed derivation [9|, are coded into the simulation 


to produce the field values for the tracks. Eq. (24) is explicitly linear in the track length 
|u|<5f, a feature we mentioned earlier. 

Eqs. (^ & ^||) are incorporated into a track-by-track Monte Carlo simulation. The 
code uses the exact formula Eq. (|2^) unless the conditions are very close to the Cherenkov 


angle, posing a 0/0 numerical problem, in which case Eq. (24) is used. 

Numerical summation of the electric field from all the tracks weighted by the proper 
charge automatically incorporates the features of coherence. As illustrated shortly, coher¬ 
ence produces a signal that peaks at the Cherenkov angle (0 C ), with a width away from 
the Cherenkov angle that shrinks with increasing frequency ||. The coherence features 
are sufficiently intricate that we have devoted separate subsections to the topic. 

5 Analysis of the effects of keeping the next order in the expansion of the phase shows that significant 
deviations from the Fraunhoffer result appear at distances where the Fresnel zone sets in jf2| , as discussed 
from a general point of view below ltd]. 
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Figure 16: Shower parametrization is done by assuming the longitudinal excess charge devel¬ 
opment to be a Gaussian n(z') moving along the shower axis (^-direction) with speed v. n is the 
direction of the observer, making an angle 9 with the shower axis, p is the transverse distance to 
the observer. 


5.2 Parametrization Method 


Charges of opposite sign which radiate coherently will give electric fields that cancel. It 
is the excess charge (e — p ) that determines the net field. Following Buniy and Ralston 
(BR) Q], one can then parameterize the excess (net) charge development in a shower 
and calculate the electric field due to the prescribed relativistic current. This is a flexible 
and compact approach that takes a few essential parameters from the shower and allows 
inspection of both the Fresnel and Fraunhoffer zones. We outline this method to keep our 
presentation self-contained. 

The vector potential A in the Lorentz gauge (adapted to e(w) and p = 1) is given by 

AMx) = —-= / d 3 x' ——— / d t'e lut J(t'x') (25) 

c\J2tt J \x — x'\ J 

The BR method parametrizes the excess charge development in the shower by a current 
along the shower axis (z) as 

J(t',x') = vn(z')f(z' - vt,\pt), (26) 


where p = ^Jx 2 + y 2 is the radial distance from the shower axis, and n{z !) is the excess 


charge (e — p ) distribution, approximated by a Gaussian near the shower maximum (see 

(27) 


Fig. jl6|), as 


n ( z ) = !h^ e - 


V2^ 


This is analogous to Greisen and Rossi’s parameterization |44j], where ?i max is the excess 
charge at the shower maximum, except that the excess (not total) charge (e—p) is modeled 
in the BR case. The parameter a corresponds to the “longitudinal spread” of the shower 
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near maximum and can be extracted by fitting a Gaussian to simulations of the excess 
charge profile of a shower. 

The region of the shower over which the fields arrive nearly in phase at a distance R(t) 


is called the coherence zone [43], and is given by 


^~coh — 


R 


k sim 8 


2 /)’ 


(28) 


where k = 2ir/\ is the wave number and 8 is the angle between R and the shower axis 
( z ). The value of the dimensionless ratio 77 , 


V = 




(29) 


determines how one should calculate the field. The ranges 77 < 1 and 77 > 1 correspond to 
the Fraunhoffer and Fresnel zones, respectively. In the case of 77 > 1, the series expansion 
of the phase, Eq. ( p0| ) fails, and the exact phase must be kept. As noted earlier, both the 
Fraunhoffer and Fresnel zones are far-field problems in the sense that kR 1 is assumed 
in both cases. 

For lo ~ GHz the Fraunhoffer approximation^ is appropriate for R > 100 m distances 
typical of the RICE experiment. However, for a test beam experiment, where the interplay 
between different length scales is important, the Buniy-Ralston method to calculate the 
field in the Fresnel zone is required. 


5.3 Radiated Power and The Form Factor 

We return to the ansatz Eq. (|2^), which describes a current with transverse distribution 
independent of the longitudinal evolution. This is a reasonable assumption after the first 
few radiation lengths of the shower development. We make the specific ansatz 

= e z qv n(z') g{f>) S(z' — vt'). (30) 


Note that the transverse extent of showers is small compared to the longitudinal extent. 
The phase and the denominator in Eq. (|25|) can be safely expanded in cylindrical coordi¬ 
nates as 


\x — a?| = \J{z — z') 2 + p 2 — 2p ■ p 1 + p' 2 

» fl(z ’*' ) ( 1_ 7F 00s *') +O (f?r) 

and 

1 _ 1 
|x — x'\ R 

6 Shown in detail earlier (Eqs. (^3] & |i^)). 


(31) 


(32) 
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where R(z, z!) = \J{z — z')' 2 + p 2 and R' = y/z' 2 + p' 2 . This expansion of the transverse 
variables avoids the Fraunhoffer approximation, and leads to a “factored” expression for 
the field: 

RE^^x) « I F(oj). (33) 


where I is an integral characteristic of the longitudinal shower history, and F(uj) is called 
the form factor. The integral I can be evaluated in the saddle-point approximation, as 
detailed elsewhere [43], or evaluated by Monte Carlo. This feature of factorization to 
express the field using a form factor works in both the Fraunhoffer and Fresnel zones. 

In the Fraunhoffer zone, there is a further, and very remarkable, simplification. Keep¬ 
ing the term linear in z\ and dropping terms of order R' 2 /R 2 in the phase expansion, the 
electric field at the Cherenkov angle (0 C ) is 


RE w (x) = q2y/2iri 


te 


,ikR 


LO 

a 


k cos 


k sin 0 r 


- e n 


x / dz'n(z') / dp'p'g(p') / d<ff e 


ev / ev 

fp—i p 1 sin 6 »c cos <t>’ 


(34) 


It is remarkable that all terms in the phase exponent depending on z 1 have vanished. This 
means that radiation from the shower is coherent over the entire length of the shower 
in the Fraunhoffer approximation: the “coherence zone” is limited only by the track 
length])]. In obtaining these results, the Cherenkov condition cos0 c = c/nv = 1 /n(3 and 
the approximations: cos0 c = z/R, sin 9 C = p/R. (see Fig. [l6|) and p'/R 1 were used. 

Inspecting the transverse integral, we see that it alone contributes to the form factor 
F(u), which is given by 

F(u) = J d 2 x' ± e- i ^<g(x l ± ), 

= / d p'p'g(p') d ^ e -^p'sm 0 ccos^ (35) 

Jo Jo 

= 2vr J dp'p'g(p') J Q (^-p'sinO^j . 


Here J D is the Bessel function of order unity. Just as in particle physics usage, F(co) is the 
Fourier transform of the (transverse) excess charge distribution. 

The preceding analysis leads to a convenient expression for the electric field at the 
Cherenkov angle in the Fraunhoffer approximation, which apart from some normalization 
factors is given by 

\RE w (x)\ = q 2a/27 rl(a) UJS ^^ C |K(a;)|. (36) 

c 

The radiated power at the Cherenkov angle (0 C ) depends on the form factor only. As noted 
earlier, the same form factor can be used in both the Fraunhoffer and Fresnel limits, so 
that laboratory information about the form factor can be used directly. The form factor 
is the central topic of subsequent sections discussing the coherence. 

'The importance of the coherence zone is well known in the Landau-Pomeranchuk-Migdal 
effect, where the lengthening coherence zone suppresses bremsstrahlung. The Cherenkov condition causes 
the coherence zone to expand to equal the entire track length. 
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6 Electric Field Pulse Calculation 


In this section we carry out the computational outline just developed. The vector nature 
of superposition is taken into account. The analytic construction of the form factor, de¬ 
veloped earlier, is adapted to the numerical calculation, so that two independent methods 
can be compared. We observe using the simulation data that coherence extends far above 
the frequency regime anticipated from simple estimates using characteristic scales such as 
the Moliere radius. 


6.1 Vector Superposition 


To calculate the electric pulse from a GEANT generated shower, we summed the con¬ 
tributions to the electric field in Eqs. (23 & ^4|) from all charged track segments us¬ 
ing full 3-dimensional geometry. We denote the observation point by a unit vector, 
h = (sin 0 cos 0, sin0sin(/>, cos0). The starting time t\ of each track segment is ob¬ 
tained from GEANT output. The time interval 5t for each track-segment is calculated for 
a particle of total energy E (from GEANT output) traveling at constant speed j3 through¬ 
out the track segment. Slowing of low energy particles is thus taken care of approximately. 
The electrical held has an azimuthal symmetry which will be investigated later. 

The electric held is a complex vector quantity. As such there are two questions of 
coherence, one to do with phase and one to do with the vector nature of the held. The 
vector nature of the held depends on the term n x (h x v) = —v± in Eqs. (23 & [24] j. 
This construction picks out the component of the velocity which is perpendicular to the 
direction of the observer. Since the track segments vary in direction, we add the electric 
held contributions as vectors which allow any cancellation that may occur. The electric 
held amplitude is then proportional to the track length transverse to the direction to 
the observer (|u|_|_(It). On the Cherenkov cone, one then expects that the held generated 
by the entire shower is proportional to the track length projected along the shower axis 
times the sine of the Cherenkov angle. The projected track length thus accounts for the 
vector nature of the electric held for tracks that go in different directions, and serves as 
an important diagnostic of the whole procedure. 

However, the projected track length says nothing in itself about the conditions for 
phase coherence, which is the topic of the next subsection. 


6.2 Phase Coherence 

The numerically generated phases of the helds, track-by-track, are determined by the 
complex exponentials in Eq. (]23|). We call the phase angles oj(t\ — • r \) and co5t(l — 

h ■ (3n) the translational phase (TP) and the Cherenkov phase (CP) respectively. The 
translational phase is kinematic, a consequence of translational invariance, and depends 
on the beginning of each track segment. The Cherenkov phase vanishes at the Cherenkov 
angle, which a point of stationary phase and dominates the emission. Two or more track 
segments will contribute in phase if the observer lies on the Cherenkov cone, and the 
beginnings of the tracks do not destructively interfere. 
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Figure 17: Distributions of the translational phase (TP) and the Cherenkov phase (CP) at the 
Cherenkov angle (6 = 6 C ) for a single shower of energy 100 GeV and threshold of 0.611 MeV. 
The TP (left plot) shows strong coherence (sharp peak) at low frequencies. A slight positive 
enhancement of the TP can be interpreted as particle positions slightly lagging the light cone 
because (3 < 1. The CP (right plot), in comparison, remains strongly peaked at all frequencies, 
which underscores the argument that the phase (CP) varies a little over the whole track. 


We studied the distribution of the TP and CP at the Cherenkov angle (9 « 55.8°) and 
at an angle (9 = 40°) off the Cherenkov angle ( Figs. 0 & [18]). We used a single 100 GeV 
shower with 0.611 MeV total energy threshold to make the phase distribution plots. We 
studied the distribution at frequencies: 10 GHz, 5 GHz, 1 GHz and 500 MHz. 

The distribution of the TP (Fig. 0 ) shows a strong peak at low frequency, indicating 
coherent phase emission from track segments at the Cherenkov angle. At high frequencies, 
the distribution tends to become flat, with random phases indicating a loss of coherence. 
A net positive phase indicates the particle positions slightly lagging behind the light cone 
due to their speed (3 slightly less than 1. The flat distribution of the TP in Fig. 18, on 


the other hand, clearly indicates random phases coming from track segments at angles off 
the Cherenkov cone. 

The distribution of the CP remains qualitatively the same both at the Cherenkov 
angle (Fig. 0 ) and at angles off the Cherenkov angle (Fig. ^). It is also apparent that 
frequency dependence of the CP is very weak. 

The preceding discussion makes clear that phase coherence is dominated by the TP. 
Finally, the distribution of TP and CP variables are substantially uncorrelated (see Figs. 
[l9| &; ^). These features indicate a valid factorization of the electric field, leading to an 
independent motivation for the Ansatz Eq. (|3^) and consequent appearance of the form 
factor. We discuss this factorization in the following subsection. 
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Figure 18: Distributions of the TP and the CP at an angle (9 = 40°) off the Cherenkov cone for 
a single shower of energy 100 GeV with total energy threshold of 0.611 MeV. The flat distributions 
of the TP (left plot) at all frequencies indicates the randomness of phases coming from randomly 
located track segments. The CP distribution (right plot), is qualitatively the same as the case on 
Cherenkov angle. Taken together, these features indicate that phase coherence is dominated by 
the TP. 
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Figure 19: Scatter plots of the translational phase (TP) and the Cherenkov phase (CP) at the 
Cherenkov angles (0 C ) and at frequencies: 1 GHz (left plot) and 5 GHz (right plot). A correlation 
between the two variables would appear as a slanted line, or similar feature. The plots demonstrate 
that the translational and Cherenkov phases are substantially uncorrelated. 
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Figure 20: Same as Fig. [l9]at angle (9 = 40°). 


6.3 Factorization of the Electric Field, and the Discrete Form Factor 


The electric field from all the shower particles can be factorized because the TP and CP 
phases are uncorrelated, as shown above. We write the field in factored form as 

TTitOt V'' *0T P / -I \ V'' l<t> TP V“~V i(jF P -IN /nry\ 

K oc 2^ e (e - 1) ~ ^ e 3 zL( e 3 - X ) . ( 37 ) 

3 3 3 


where the sums are over all track segments. The factorization takes a rather simple form 
at the Cherenkov angle. The electric field close to the Cherenkov angle in Eq. (24) can 
be written approximately for a single track as 


RE U (x) = (±\ v±6t e ikR e ^(i-fcosd) e -Offi-r ix. ( 38 ) 

V2vr \c z J 

The approximation t\ ~ z\/v has been assumed in writing Eq. (|38|), which is well sup¬ 
ported by the TP phase plot at the Cherenkov angle shown above. The phase angle 
loz\{1/v — ncos 9/c) vanishes for the Cherenkov condition (cos0 c = c/nv) and we have full 
coherence along the ^-direction (shower axis) as observed earlier. 

We now write the factorization in Eq. (^) for the total electric field at the Cherenkov 
angle from all the charged particles in the shower as 


RET{x) 


iu; ( q \ 

V2tt WJ 


e ikR Y J S j {vi_) j 5t j e _i2 w "'(riuP 
3 


(39) 


where s- 7 = ±1 for positrons and electrons respectively. The total electric field in Eq. (|39| ) 
is thus proportional to the total track length (fffj s 3 [y±) 3 St J ) transverse to the direction 
of the observer at any frequency lo. This track length is the projected (e — p) track length 
times sin 6 C , as described earlier. 

The coherent electric field emission at different frequencies now depends on the Monte 
Carlo “discrete form factor” F(uj) mc, given by 

F(u) MC = Y / s j e- iS ? fl -^ j =Y J s j e~ i ^< M . (40) 

3 3 
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Figure 21: Absolute value of the form factor (|F(u;)|) (left plot) and the frequency spectrum 
(w|i 7 '(w)|) (right plot) plots for 1 TeV, 500 GeV and 100 GeV showers. The form factor has been 
calculated using Eq. ([h]), where we naively sum over all the shower particles along the z-axis 
(shower axis). The curves are normalized by dividing by the particle count N. 


Note that the observation point is in the x — z plane. 


6.4 The Discrete Form Factor and the Frequency Spectrum 


We calculated the discrete form factor in Eq. ( |40| ) for 1 TeV, 500 GeV and 100 GeV 
showers. A factor of lo times the form factor is proportional to the electric field amplitude. 

First, we naively carried out the sum in Eq. (^0|) for all the shower particles along 
the whole shower axis ( z ) as suggested by the Cherenkov condition in the Eq. (|38|). The 
absolute value of the form factor and the frequency spectrum are plotted in Fig. 21. It 
shows an extended region of coherence, with \Formfactor\ ~ N/uj, where N is the number 
of shower particles. This naive estimate is somewhat misleading because the Cherenkov 
condition is not satisfied for the whole shower axis. The CP also becomes flat at very high 
frequencies, as does the TP at the Cherenkov angle (see Fig. 0a)- A further complication 
is statistical fluctuations, which become large where the form factor is small. 

For a better understanding of the form factor and the frequency spectrum, we now 
turn to the analytic method. 


6.5 The Analytic Form Factor and the Frequency Spectrum 

The analytic form factor as defined in Eq. (^) is the Fourier transform of the snapshot 
of the charge distribution. It is a good approximation to assume that most of the electric 
field contribution comes from the particles at the shower maximum [|43[ j. 

To calculate the analytic form factor, we then determined the transverse distribution 
of the particles within half a radiation length on both sides of the shower maximum for a 
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Figure 22: Radial distributions dN/dp of the excess (e — p ) electrons within half a radiation 
length on both sides of the shower maximum for a 100 GeV shower (averaged over 50 showers). 
The solid curve is a Pade(l,3) fit to the excess charge distribution. 


100 GeV shower (averaged over 50 showers). We found both the distribution of the total 
particles d N(e + p)/dp and the distribution of the excess electrons d N(e — p)/dp. Here, 
p is the cylindrical radial distance from the shower axis. The excess charge distribution 
(dJV(e — p)/dp) is plotted in Fig. |22. 

We observed that the radial distribution falls approximately like 1/p 2 for large p, 
and peaks at a remarkably small value of 0.5 cm. One way (by no means complete) to 
understand these results recalls the multiple scattering of an electron by a static Coulomb 
held. For a screened Coulomb held, the modihed Rutherford scattering cross section is 
given by [39]: 


dus 

d Q 


= 87T 


Ze 2 


Q 


pc j ( Q 2 + Qir 


(41) 


where Q is the momentum transfer and Z is the atomic number of the scatterer; Q s = 
(Z 1 / 3 /192)mc is the momentum transfer associated with the screening radius; m is the 
mass of the electron. For elastic scattering in the ultra-relativistic limit, Q 2 = 2p 2 (l — 
cos 9) ~ p 2 9 2 , where p is the electron momentum. The cross-section formula in Eq. © 
now takes the form: 


cK / Ze 2 \ 2 p 2 6 

d 0 ^ \ Pc ) ( p 2 0 2 + Q2)2- 


(42) 


The average deflection angle, <#>« p/Az. We can get a crude estimate of the 
transverse distribution dN/dp from Eq. (|4^) as: 


dlV / Ze 2 \ 2 (nAz 3 \ p 

dp ~ ^ 1 k (3c) ^ p 2 J (p2 _|_ p 2f 


(43) 
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where n is the number density related to N by the formula N = noAz. The peak of this 
distribution is given by 

ft = ^4 (44) 

P 

We evaluate this with the average energy of the particles taken to be approximately 
equal to the critical energy (E c ) at the shower maximum. The peak of the distribution is 
then p Q = 0.5 cm for ice (Z = 7.2, A z = 39 cm and p ~ E c = 70 MeV). The excellent 
agreement may be fortuitous, but gives some confidence that the Monte Carlo excess 
charge distribution, which also includes numerous atomic collision, pair production and 
Compton scattering processes, has a simple physical origin. 

We made a fit to the excess charge distribution with a Pade (1,3) approximation of the 


form dN/dp = f(p) = n(p — o)/(l + bp + cp 2 + dp 3 ) (see Fig. 22). The fitting parameters 
are: n = 2.15, a = 0.01, b = 4.73, c = 6.96, d = 0.33. The fit is good up to p ~ 20 cm. 
The choice of the Pade fit was made to preserve the p 1 geometric zero in dN/dp at the 
origin, and the 1 /p 2 asymptotic behavior. 

Finally, we calculated the analytic form factor using Eq. 


as: 


F(uj) = 2n J o fip) J o P sin 0 c ^j ■ 


(45) 


For the high frequency end it was necessary to use a convergence procedure to modulate 
the Bessel transform. We compare the analytic prediction to the discrete form factor of Eq. 
(|4C|) at the shower maximum for 100 GeV, 500 GeV and 1 TeV showers. The comparison 
is shown in Fig. [2^. For reference we also plotted the analytic and the discrete frequency 
spectra, which are lo times the absolute value of the form factor. The agreement between 
the two methods is good up to 10 GHz for all energies and up to 50 GHz for 1 TeV. 


6.6 Direct Calculation: Monte Carlo Field Spectrum 

The frequency spectrum of the electric field calculated using Eqs. (|23] & |1|) at the 
Cherenkov angle and an angle off the Cherenkov angle are plotted in Fig. |24|a. We 
calculated the spectrum for 1 TeV, 500 GeV and 100 GeV showers (each averaged over 
many showers). The electric field amplitude rises linearly at low frequency: this is the lin¬ 
ear dependence on ui due to dimensional analysis. The coherent behavior at the Cherenkov 
angle and the incoherent behavior at a widely separated angle (40°) are clear. The fre¬ 
quency spectrum of the electric field calculated from a 100 GeV GEANT generated shower 


is compared to the same from the ZHS code in Fig. 24 b. 

From the previous analysis, the frequency dependence of the electric field amplitude is 
expected to be u>\F(u)\. We compare u;|-F(u;)| to \E U \ generated by the Monte Carlo in Fig. 
|25| . The electric field amplitude at the Cherenkov angle is proportional to the projected 
(e — p) track length x sin 0 C as shown in Eq. (fkil). We used the total projected (e — p) track 
length value of 70 meter (as we found before) for a 100 GeV shower as a pre-factor to 
normalize the analytic form factor. It is clear from the plot that the analytic form factor 
and simple dimensional considerations explain the electric field amplitude rather well up 
through 5-15 GHz range. The coherence persists up to frequencies as high as 50 GHz. We 
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Figure 23: Absolute value of the form factor |F(u;)| (left) and the frequency spectrum w|F(w)| 
(right) plotted versus frequency for 1 TeV, 500 GeV and 100 GeV showers. The discrete form 
factor has been calculated using Eq. (ff0|), where we sum over all the particles within a radiation 
length at the shower maximum. The analytic form factor has been calculated from the Fourier 
transform to a fit to the excess charge distribution (see Fig. |2^). Both curves are normalized by 
dividing by the particle count N Note that the analytic curve tends to be lower than the numerical 
calculation at high frequency. 


notice that the agreement between the direct Monte Carlo calculation of the frequency 
spectrum and the analytic frequency spectrum gets better with increasing shower energy. 
We believe that the analysis captures the important physics of the processes, and validates 
the results of the Monte Carlo. 

6.7 Related Issues 

Electrons and positrons undergo Coulomb scattering while traversing through the medium. 
The particle track is therefore not a straight line from the point of its creation to the point 
where it falls below threshold. The track contains many kinks due to elastic Coulomb 
scattering. GEANT gives a detailed output of particle tracks which contain these kinks as 
stated earlier; we calculate the electric field from step-tracks. The Monte Carlo developed 
by ZHS, on the other hand used the straight tracks from the start points to the end points, 
with timing correction, to calculate field Q to simplify their calculation. 

We studied the effect of these kinks and also the effect of taking all tracks along the 
shower axis (see Fig. |26| ). We found little difference between the cases when we included 
(the usual case) and not included the kinks. This is a result of the extended coherence 
zone along the shower axis at the Cherenkov angle. When we calculated the electric field 
taking only the components parallel to the shower axis, we found that the field increases 
almost linearly with frequency, as expected from a single charge. 

For practical purposes, particles are removed from the Monte Carlo simulation once 
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Figure 24: Direct Monte Carlo calculation of the frequency spectrum of the electric field mag¬ 
nitude for 1 TeV, 500 GeV and 100 GeV showers (left). The calculation is done using Eqs. (£3] 
& HH) at the Cherenkov angle and at 40° angle. The coherent behavior at the Cherenkov angle 
and the incoherent behavior at the 40° viewing angle are clear. The frequency spectrum at the 
Cherenkov angle for a 100 GeV shower is also compared to the same calculated using ZHS code 
(right). All thresholds are 0.611 MeV and the averages are done with 50 showers in case of 100 
GeV, 20 showers in case of 500 GeV and 10 showers in case of 1 TeV. 



Figure 25: Comparison between the frequency spectra of the electric field magnitude at the 
Cherenkov angle from direct Monte Carlo calculation and from analytic calculation. The analytic 
calculation is done using Eq. (40) with the total projected track length of 70 m for a 100 GeV 
shower as found before. 
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Figure 26: Frequency spectrum of the electric field magnitude at the Cherenkov angle for a 
100 GeV shower (averaged over 50 showers). Our study shows that the electric field calculated 
using particle tracks with and without kinks make little difference. The electric field magnitude, 
on the other hand increases almost linearly with frequency like a single charge when we took the 
component of the tracks parallel to the shower axis only. 


they fall below certain threshold energy (0.611 MeV in our case). While physical particles 
do not suddenly stop moving, simulated particles in the Monte Carlo, have potential 
emission of stopping radiation. Even though high energy particles are almost parallel to 
the shower axis, the low energy particles close to threshold are isotropic. The stopping 
radiation (which is forward) from all these isotropic particles in general could cause a 
problem if they are aligned with the observer on the Cherenkov cone. 

Our study shows (see Fig. 27) that there is indeed an isotropic component in the 


shower below energy 0.85 MeV. However, they do not contribute to the radiation at the 
Cherenkov angle ( 9 C ~ 55.8°). Particles with energy 1 MeV, for example, Cherenkov 
radiates at an angle 49°. 


6.8 Angular Pulse Distribution 

The angular distribution of the electric field amplitude peaks at the Cherenkov angle 
(0 C ~ 55.8°). Fig. [2^ shows the 1 GHz pulse from a 100 GeV shower (averaged over 
50 showers) with 0.611 MeV total energy threshold at 4 different azimuthal angles cj) = 
90°, 180°, 270° and 360°. The pulse shows very little dependence on cj) which corresponds 
to the approximately symmetric distribution of particles about the shower axis as stated 
earlier. We have also plotted the electric field amplitude from the ZHS Monte Carlo for a 
100 GeV shower (averaged over 50 showers) at 1 GHz frequency in Fig. ^ for comparison. 

Fig. |2^ shows the angular pulse distribution for a 100 GeV shower (averaged over 
50 showers) at 1 GHz, 750 MHz, 500 MHz and 250 MHz frequencies. The Gaussian half 
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Figure 27: Plot of the cosine of the polar angle (0) of each track-segment with respect to the 
shower-axis versus its total energy. Although high energy particles are very much aligned with 
the shower axis, there is an isotropic low energy component which moves with the shower. The 
isotropic component of the shower does not contribute to the Cherenkov pulse, which validates the 
removal of the particles when they fall below a certain low-energy threshold, say 0.611 MeV. 
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Figure 28: Angular pulse distribution of a 100 GeV shower (averaged over 50 showers) with 0.611 
MeV total energy threshold from GEANT. The pulse is calculated at 1 GHz frequency and at 4 
different azimuthal angles (c f >). The Cherenkov peak at the observation angle 0 = 0 C shows very 
little dependence on </>. The same pulse from the ZHS Monte Carlo is also plotted for comparison. 
















Figure 29: Angular pulse distribution of a 100 GeV shower (averaged over 50 showers) with 
0.611 MeV total energy threshold from GEANT. The pulse is calculated at frequencies: 1 GHz, 
750 MHz, 500 MHz and 250 MHz. The plot shows inverse scaling of Gaussian half width of the 
Cherenkov pulse with frequency. 


width of the pulses are 2.0°, 3.7° and 7.3° at 1 GHz, 500 MHz and 250 MHz frequencies re¬ 
spectively. This corresponds to approximate inverse scaling of pulse width with frequency 
which is analogous to a single-slit diffraction pattern as pointed out in reference pj]. 


7 Summary of Results and Conclusions 


We have analyzed 100 GeV - 1 TeV electromagnetic showers and the radio frequency radi¬ 
ation they produce in great detail. These studies are a necessary ingredient for designing 
an experiment for radio detection of UHE cosmic ray induced showers in radio- transpar¬ 
ent media. In particular, experiments to detect radio emission from showers induced by 
high energy cosmic ray neutrinos interacting in the surface of the moon [ 451 and in the 
South Polar ice-cap |4(| are underway and have reported preliminary results. Coherent 
radio emission from electromagnetic showers has recently been demonstrated in the labo¬ 
ratory (T|]. The technique is gaining recognition as a powerful tool for particle detection. 
The thorough dissection and understanding of all the intricacies of the showers and their 
relationships to the final radio pulse produced is our goal achieved in this paper. 

Energy information for each stage of every track is readily available from the GEANT 
shower code. This information is essential for tracking the energy distribution in the 
shower and identifying the sources of radio emission. Among our results, is the direct 
determination of the radiation length in ice from an exponential fit to bremsstrahlung 
radiation energy loss as a function of depth in Sec. 3.1. When this information is combined 
with the direct extraction of the ionization loss, Sec. 3.3, we determined the critical energy 
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from the data shown in Fig. |j. The value obtained is nicely consistent with that found by 
using the Moliere radius extracted from the data for the radial energy flow, Sec. 3.2, in 
combination with the radiation length. The consistency between the critical energy values 
indicates that our application of the code and analysis of the data gives a correct physical 
picture of the interplay among the competing processes in the shower as it develops. 

Because the coherent radio emission of interest depends upon the charge excess in the 
shower, we need the energy profile of the contributions to the charge excess. This again 
requires the GEANT track-by-track energy information, and we show the total charge 
imbalance broken down into energy ranges in Fig. 12. A large fraction, more than 50%, of 


the imbalance comes from the energy range below 5 MeV. Though the track lengths are 
small, the large number of particles leads to a significant contribution to the total track 
length of the shower. 

The total and projected track lengths for both the total and excess charge populations 
are determined and shown to be proportional to total shower energy in Sec. 3.5. The total 
track length is comfortably below a rough estimate of the upper bound for the total track 
length. Our detailed study of the longitudinal profiles in iron shows overall agreement but 
differences in detail among the GEANT, the EGS4 and the ZHS simulations in Fig. |8[ 
The profiles in ice for 100 GeV, 500 GeV and 1 TeV electron and photon induced showers 
in ice are well described by a modified Greisen parametrization with critical energy value 
extracted from the simulation data, as described above, and two fit parameters. The 
confidence levels of the fits are typically 80% - 90%, as summarized in Tables 1 Sz 2. 
The GEANT profiles in ice are qualitatively similar to those from the ZHS code, but 
lie typically 25% - 35% lower. To determine whether small differences in cross section 
values used in the different simulations could account for this difference, we developed a 
1 -dimensional shower code with the full set of cross sections for the relevant processes. As 
we report in Sec. 3.7, the depth at maximum and the number of particles at maximum are 
rather insensitive to changes in the cross sections. We therefore believe that the differences 
in cross-sections are not the source of profile and total track length differences between our 
simulation and that of ZHS. Unfortunately, one important check that we have not been 
able to make is the GEANT vs. ZHS ionization energy loss (dE/dx). The ZHS code does 
not admit a readout of energy loss by shower particles. Therefore, we were unable to make 
from ZHS code plots similar to Figs. [3] Sz § we made with GEANT. If the ionization loss 
in ZHS code is much lower than in GEANT, it might account for the difference between 
showers produced by the two Monte Carlos. 

We develop the framework to calculate the electric held in the 100 MHz to multi 
GHz frequency range in Sec. 5. The track-by-track held calculation, applicable in the 
Fraunhoffer zone, and the calculation treating the shower as a continuous current density, 
applicable in the Fraunhoffer and and Fresnel zones, complement each other. We present 
both for this reason. The direct, numerical calculation of the total held from the vector 
sum of the helds from individual tracks allows a detailed study of the dependence on to¬ 
tal track length. Moreover, we elucidate the effects of random direction changes due to 
collisions and the pattern of phase relationships among the contributions from all tracks. 
Complementing these intensive numerical calculations, we present analytic work that em¬ 
ploys the form factor of the effective current density to calculate the held in both Fresnel 
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and Fraunhoffer regimes. Given our realistic model for the current density, which we fit 
to the transverse charge distribution from the simulation, we show the remarkable result 
that the radiation from the shower is coherent over the whole shower at the Cherenkov 
angle in the Fraunhoffer limit. 

Our study of the track-by-track phase coherence in Sec. 6 reveals that the phase 
associated with the initial time and position of a track (TP) and the phase associated 
with the “diffractive emission” (CP) from the track as a whole are uncorrelated. This 
result supports our model of the current, and consequently of the field, in Sec. 5. Figs, 
m show the strong phase peaking of the “diffractive” phase, called Cherenkov phase 
(CP) in the discussion, on and off the Cherenkov angle as defined by the shower axis. 
The phase associated with the initial coordinates of the track, called the translation phase 
(TP), is coherent at the Cherenkov angle but completely random away from that angle. 

We noticed some time ago [[47]] that the frequency spectrum of the electric field at 
the Cherenkov angle, calculated directly from the track data and the Fraunhoffer zone 
formulas developed in Sec. 5., flattens out at frequencies above 2 GHz, as shown in Fig. 
24. It is clear from the figure that the ZHS simulation shows the same behavior but 
they did not examine the high frequency region further |J. We employ the form factor, 
which correctly accounts for the transverse spread of the shower and its effect on the 
field to analyze the frequency spectrum. We evaluate its Fourier transform, F(u>), two 
independent ways and find agreement. We then show that the observed spectrum of the 
field is faithfully represented up to 5 GHz for 100 GeV showers and up to 15 GHz for 500 
GeV and 1 TeV showers. Up to 15 GHz, we now have a clear picture of the behavior of 
the frequency spectrum, which continues to show coherent behavior. 

The electric held in the Fraunhoffer zone as a function of the observation angle from 
the shower axis peaks at the Cherenkov angle. The width of the peak shrinks inversely 
with frequency and the height (held strength) rises linearly with shower energy. These 
features confirm the ZHS results and the off-Cherenkov angular dependence is also similar. 
The height of the peak at the Cherenkov angle for a given energy and frequency is larger 
in their case by about 25-35%. This is perhaps not a surprising difference between two 
independent shower simulations and held calculations. We have performed an extensive 
set of tests and cross-checks to validate our results. 

7.1 Conclusions 

Our study quantihes coherent Cherenkov radiation from high energy showers and shows 
that coherence persists from 100 MHz to tens of GHz. The existence of Coherent Cherenkov 
radiation goes back to Askaryan, and has been studied for decades, while the persistence 
of coherent emission at the Cherenkov angle at multi-GHz frequencies is new. The persis¬ 
tence of coherence is established by our new, highly detailed study of the actual shower 
currents and corresponding phase distributions in CP and TP introduced in Sec. 6. The 
recognition that the TP is coherent over the whole length of the shower at the Cherenkov 
angle is a new insight into the connection between shower particles and the fields they 
produce. We have analyzed the electromagnetic shower characteristics in great detail, 
including the energy distributions and energy losses, so that a complete set of shower 
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parameters can be extracted from the simulation data to confirm and cross check the con¬ 
sistency of our picture. A grand summary of pertinent parameters and comparisons with 
the ZHS and other sources where available is presented in Table 5. 


Table 5: Comparisons between 100 GeV showers in ice from GEANT 
and from ZHS Monte Carlos 


Parameter 

GEANT 

ZHS 

Total Energy Threshold (MeV) 

0.611 

0.611 

Total absolute track length (meter) 

399 ± 5 

642 

Total projected (e + p) track length (meter) 

374 ± 4 

519 

Total projected (e — p) track length (meter) 

70 ± 8 

131 

Position of the shower max. (radiation length) 

6.5 

7.0 

Number of particles (e + p) at shower max. 

111 ± 7 

155 ± 10 

Excess electrons (e — p) at shower max. 

20 ± 2 

37 ± 3 

Fractional charge excess at the shower max 

~ 18% 

~ 24% 

Cherenkov peak at 1 GHz (Volts/MHz) 

7.5 x 10 -9 

1.0 x 10 -8 


The simulation and field calculation developed are directly applicable to energies up 
to 1 TeV. Higher energy showers can be ’’bootstrapped” by evolving the shower particles 
into this regime with a corresponding multiplication in particle number and emitted field 
intensity. We intend to parametrize and extrapolate the data to the multi-PeV data 
needed for the UHE showers relevant to the RICE experiment. With the GEANT base, it 
will also be interesting to investigate similar shower and field calculation for the hadronic 
component of neutrino-induced showers. 
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